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Abstract 

A  three-dimensional  (3D)  numerical  model  associating  the  heat/mass  transfer  and  the  electrochemical  reaction  in  a  proton  exchange 
membrane  (PEM)  fuel  cell  is  developed  in  this  study,  and  a  miniaturized  PEM  fuel  cell  with  complex  flow  channels  is  simulated.  The 
numerical  computation  is  based  on  the  finite-volume  method.  Governing  equations  for  flow  and  heat/mass  transfer  are  coupled  with  the 
electrochemical  reactions  and  are  solved  simultaneously.  The  latent  heat  from  the  condensation  of  water  vapor  in  cathode  channel,  if  any,  is 
considered.  The  perimeters  of  the  bipolar  plates  are  also  included  in  the  computational  domain  to  account  for  their  heat  conduction  effect.  The 
miniaturized  PEM  fuel  cell  has  a  membrane  electrode  assembly  (MEA)  sandwiched  by  two  brass  bipolar  plates  etched  with  a  number  of 
winding  gas  channels  with  a  flow  area  of  250  pm  x  250  pm.  The  influence  of  anode  gas  humidity  on  the  performance  of  the  fuel  cell  is 
investigated  through  model  prediction.  Finally,  field  details  of  velocity,  mass  fraction  and  electromotive  force  are  illustrated  and  discussed. 
©  2003  Published  by  Elsevier  Science  B.V. 

Keywords:  PEM  fuel  cell;  3D  numerical  model;  Complex  flow  channels;  Heat/mass  transfer 


1.  Introduction 

Proton  exchange  membrane  (PEM)  fuel  cells  operate  at  a 
low  temperature  (from  room  temperature  to  90  °C)  to  con¬ 
vert  chemical  energy  directly  into  electrical  energy.  Low 
emission  of  pollutants  and  high  electricity  conversion  effi¬ 
ciency  make  this  technology  attractive,  and  it  is  currently 
being  considered  as  one  of  the  most  promising  candidates 
for  stationary  and  automotive  applications  as  a  substitute  of 
traditional  power  systems  (i.e.  thermoelectric  power  plants 
and  internal  combustion  engines)  [1-4]. 

A  PEM  fuel  cell  consists  of  three  major  components:  an 
anode,  typically  featuring  a  platinum  or  platinum-containing 
catalyst;  a  thin,  solid  polymeric  membrane  which  acts  as 
the  electrolyte;  and  a  cathode,  also  platinum-catalyzed.  The 
reactions  in  a  PEM  fuel  cell  fed  with  a  hydrogen-containing 
anode  gas  (hydrogen  with  some  amount  of  water  vapor)  and 
an  oxygen-containing  cathode  gas  (air)  are: 

EL  2H+  +  2e_  (anode)  (1) 

2H+  +  2e“  +  \  02  H20  (cathode)  (2) 

H2+i02~H20  (overall)  (3) 
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From  Nernst  equation,  the  thereby  induced  electromotive 
force  or  electrical  potential  between  the  cathode  and  anode 
is  a  function  of  the  temperature  at  the  reaction  site  and  the 
partial  pressures  of  reactants  and  products  as  well: 


~  AG°  RT  In  ( Pu^  \ 
2 F  2 F  (  PHlo  ) 


where  AG0  is  the  standard  Gibbs  free  energy  change  (high 
heat  value)  of  the  electrochemical  reaction  at  temperature  T\ 
the  P  with  the  respective  subscripts  are  the  partial  pressures 
of  the  corresponding  reactants  and  products  at  the  electro¬ 
lyte/electrode  interfaces;  F  the  Faraday’s  constant;  and  R  is 
the  universal  gas  constant. 

As  typical  component  of  a  PEM  fuel  cell,  a  plate  (made  of 
high  conductivity  material)  etched  with  winding  grooves 
serving  as  gas  channels,  is  shown  in  Fig.  E  Two  such  plates 
sandwiching  a  layer  of  membrane  electrode  assembly 
(MEA)  on  their  grooved  side  can  compose  a  single  PEM 
fuel  cell.  Fuel  is  fed  to  the  channels  of  the  anode  side,  and  in 
a  similar  way,  air  is  fed  to  the  channels  of  the  cathode  side. 
With  many  such  single  fuel  cell  units  being  laminated  in 
series,  a  cell  stack  can  be  constructed.  In  practical  terms, 
plates  with  channels  etched  on  two  sides  are  used  to  sand¬ 
wich  the  MEA  so  as  to  construct  a  PEM  fuel  cell  stack.  Such 
a  plate  is  referred  to  as  a  bipolar  plate  [1] .  A  detailed  study  of 
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Nomenclature 
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AA 

CP 
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AG0 
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^drag 


P 
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R 
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^cell 
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x,  y,  z 

X 

Y 


active  area  of  electrochemical  reaction  (m  ) 
unit  area  (m  ) 

specific  heat  at  constant  pressure  (J/(kg  K)) 
mass  diffusion  coefficient  (m  /s) 
electromotive  force  (EMF)  (V) 

Faraday’s  constant  (96,496.7  C/mol) 
standard  state  Gibbs  free  energy  change  for 
water,  formed  from  oxygen  and  hydrogen  (J/kg) 
height  of  computation  domain  (m) 
current  density  of  fuel  cell  (mA/cm  ) 
current  (A) 

thermal  conductivity  (W/(m  K)) 

mass  flux  (kg/(m  s)) 

formula  weight  (g/mol) 

equivalent  weight  of  membrane  (kg/mol) 

electro-osmotic  coefficient 

pressure  (Pa) 

heat  flux  (W/m2) 

heat  source  (W/m  ) 

universal  gas  constant  (8.314  J/(mol  K)) 

source  term  in  equations 

temperature  (K) 

velocity  in  v-direction  (m/s) 

velocity  in  y-direction  (m/s) 

fuel  cell  voltage  (V) 

velocity  in  ^-direction  (m/s) 

Cartesian  coordinates  (m) 
mole  fraction 
mass  fraction 


Greek  letters 
P  permeability 

5  thickness  of  electrode  and  electrolyte  (m) 

s  porosity 

f]  activation  polarization  (V) 

X  water  content 

jj  dynamic  viscosity  (Pa  s) 

p  density  (kg/m  ) 

Q 

Pdry  dry  density  of  membrane  (kg/m  ) 

o  conductivity  (Q  cm)-1 

Subscripts 
a  anode 

c  cathode 

H2  hydrogen 

H20  water  vapor 

in  at  inlet 

m  membrane 

02  oxygen 

Superscripts 
a  anode 

air  air 


c  cathode 

e  effective 

fuel  fuel 

H2  hydrogen 

m  membrane 

02  oxygen 

Sat  saturation 


heat/mass  transfer  for  a  cell  stack  with  all  the  single  fuel  cell 
units  subject  to  investigation  would  be  a  difficult  task. 
However,  most  of  the  single  fuel  cell  units  laminated  in  a 
stack  can  operate  under  the  same  conditions  if  the  fuel  and 
oxidant  are  provided  in  a  parallel  style.  Therefore,  a  single 
PEM  fuel  cell  constructed  using  the  plates  illustrated  in 
Fig.  1  will  be  studied  in  the  present  modeling  work.  This  will 
provide  a  useful  understanding  of  not  only  a  single  fuel  cell 
unit  but  also  a  stack  of  PEM  fuel  cells. 

The  flow  and  heat/mass  transfer  in  a  PEM  fuel  cell  have 
some  special  features.  First  of  all,  there  are  two  flow  streams, 
anodic  gas  and  cathode  gas,  to  be  considered  in  a  PEM  fuel 
cell.  With  the  processing  of  the  electrochemical  reaction,  the 
species  concentration  and  overall  mass  flow  rate  in  both  the 
streams  experiences  associated  variation.  Therefore,  analy¬ 
sis  for  the  fields  of  flow,  temperature  and  species  concen¬ 
tration  has  to  reflect  the  strong  coupling  character  between 
the  flow  streams  of  the  cathode  and  anode  [5,6].  On  the  one 
hand,  the  species  and  overall  mass  variation  in  the  two  streams 
of  the  anode  and  cathode  are  affected  by  the  electrochemical 
reaction  proceeding  at  electrode/electrolyte  interfaces.  On 
the  other  hand,  however,  the  species  concentrations  and  tem¬ 
perature  distributions  are  decisive  factors  in  the  electromotive 
force  as  expressed  by  the  Nemst  equation,  which  therefore 
affects  the  electrical  performance  of  the  fuel  cell  [7,8]. 

In  previous  studies,  most  models  focused  on  one-dimen¬ 
sional  (ID)  and  two-dimensional  (2D)  flows  with  transport 
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Fig.  1.  Copper  plate  etched  with  winding  gas  channels. 
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of  the  reactants  and  products  in  the  flow  channels  and  across 
the  membrane.  In  particular,  those  models  only  treated 
simple  fuel  and  air  channel  structures.  This  allowed  the 
fundamentals  of  numerical  simulation  for  PEM  fuel  cells  to 
be  established.  For  example,  in  the  work  of  Springer  et  al. 
[9],  they  considered  a  pseudo  one-dimensional  model  in 
which  flows  in  both  fuel  and  air  channels  were  treated  as 
being  well  mixed.  Fuller  and  Newman  [10]  and  Nguyen  and 
White  [11]  developed  two-dimensional  heat  and  water 
transport  models  that  accounted  for  variation  in  temperature 
and  membrane  hydration  conditions  along  the  flow  chan¬ 
nels.  Both  of  these  models  assumed  a  well-mixed  concen¬ 
tration  in  flow  channels.  However,  a  fully  three-dimensional 
(3D)  solution  to  the  Navier-Stokes  equations  for  the  flow 
channels  is  important  to  allow  one  to  account  for  the  effect  of 
the  dimensions  of  the  flow  channel  on  the  velocity  distribu¬ 
tion.  Also,  such  a  full  solution  will  help  one  to  understand 
how  the  mass  consumption  in  the  electrochemical  reactions 
affects  the  momentum  transport  equations  [12]. 

As  a  common  practical  structure,  winding  fuel  and  air 
channels  are  adopted  in  PEM  fuel  cells  [13].  Due  to  the 
complexity  of  the  flow  channels,  a  three-dimensional  model 
is  necessary.  This  will  be  one  of  the  focuses  of  the  present 
study. 


2.  Model  development 

From  the  point  of  view  of  numerical  analysis,  the  concern 
is  to  predict  the  details  of  the  flow  and  heat/mass  transfer 
of  two  different  streams  of  fuel  and  air  separated  by  the 
membrane.  However,  the  heat/mass  transfer  are  strongly 
coupled  between  the  fuel  and  air  flows  through  the  heat 
conducting  gas  channel  plate  and  the  mass  diffusible  mem¬ 
brane.  The  two  channel  plates  and  membrane  need  to  be 
integrated  into  one  computation  domain.  Therefore,  the 
computation  domain  of  the  v-y-plane  includes  the  whole 
copper  plate  as  already  shown  in  Fig.  1.  The  flow  directions 
of  fuel  and  air  are  arranged  as  schematically  illustrated 
in  Fig.  2. 

Governing  equations  for  flow,  temperature  and  mass 
fractions  or  concentrations  in  the  PEM  fuel  cell  in  three 
dimensions  for  the  computation  domain  including  the  MEA, 


cathode  and  anode  plates  (each  in  thickness  of  0.525  mm 
and  both  of  them  have  gas  channels  etched,  as  illustrated  in 
Fig.  1)  are: 
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Due  to  the  variation  of  gas  species’  mass  fractions,  the  local 
fluid  properties  could  vary  significantly.  This  is  also  a 
reflection  of  the  strong  coupling  of  velocities,  temperature, 
and  mass  fractions,  since  the  determination  of  the  local 
fluid  properties  will  depend  on  the  local  mass  fractions  and 
temperatures. 

The  additional  source  terms  in  momentum  transport 
equations,  —pu/fix,  —pv/fiy  and  —pw/fiz  are  only  induced 
in  the  electrode  diffusion  layer.  The  purpose  of  this  is  to 
consider  the  influence  of  porous  media  to  flow  and  mass 
transfer.  It  is  assumed  that  the  permeability  in  electrode 
diffusion  layers  is  isotropic,  and  thus  there  is 

fix  fiy  fiz 


(13) 
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For  flow  and  mass  fraction  equations,  only  the  gas  diffusible 
region  is  subject  to  investigation.  Treatment  through  numer¬ 
ical  methods  [14]  will  be  implemented  to  induce  zero 
velocities  in  all  of  the  solid  regions  in  the  computation 
domain. 

In  the  energy  conservation  equation,  the  concentration- 
gradient-induced  heat  diffusion  is  very  small  and  thus  is 
neglected  [15,16].  The  mass  fraction  of  water  is  monitored  at 
local  locations.  In  case  the  calculated  partial  pressure  of 
water  vapor  is  higher  than  its  saturated  pressure,  it  is  judged 
that  a  condensation  of  water  vapor  may  happen  and  that  the 
condensation  latent  heat  will  be  considered  as  a  heat  source 
at  the  electrode  surface. 

Due  to  the  electrochemical  reaction  and  electro-osmotic 
effect,  mass  diffusion  of  water  and  hydrogen  across  the 
membrane  will  occur.  The  electrochemical  reactants  and 
products  must  diffuse  across  electrode  layers  to  reach/leave 
the  reaction  site,  that  is,  the  interface  of  the  electrolyte  and 
electrodes.  The  mass  diffusion  inside  the  porous  electrodes 
is  much  weaker  than  that  in  the  fluid  region,  so  effective 
diffusivities  for  mass  transport  within  the  electrodes  are 
adopted  by  inducing  a  porosity  coefficient  [17,18]: 


De  —  D  g3/2 
+/,m  ujw\c’ 


(14) 


2.7.  Local  current  density  and  determination  of  mass 
diffusion  fluxes  to  gas  channels 

The  proton  conduction  is  considered  as  the  electric  current. 
In  a  plate-type  fuel  cell,  the  electric  charge  passes  through 
the  electrolyte  layer  perpendicularly  and  is  collected  by 
the  channel  plates.  As  a  function  of  the  local  electromotive 
force,  local  current  density  can  be  calculated  by  the  following 
correlation: 

•  _  _ ~  *7 act)  ~  jjsell _  ✓  4-x 

“  OW *m)  +  (Sa/*a)  +  05c K)  1  ; 

where  a  is  the  electric  conductivity  of  the  electrodes  or 
proton  conductivity  of  the  electrolyte;  5  represents  the 
thickness  of  the  electrodes  and  electrolyte;  and  f]act  is  the 
potential  drop  due  to  activation  polarization.  In  the  expres¬ 
sion  of  Eq.  (4),  local  partial  pressures  of  gases  in  the  fuel  and 
air  channel  must  be  used  to  calculate  the  electromotive  force 
E.  Because  the  local  partial  pressures  are  obtained  from  the 
solution  of  the  governing  equations,  the  convection  concen¬ 
tration  is  thus  implicitly  considered  in  E . 

By  performing  a  summation  or  integration  of  the  local 
current,  the  total  current  output  from  the  fuel  cell  can  be 
obtained: 

I=J2(iAA)  (16) 

To  start  a  numerical  simulation  for  the  fuel  cell,  the  involved 
total  current  must  be  a  prescribed  electrical  parameter,  while 
the  cell  output  voltage  Vcen  and  the  local  current  density  will 
be  calculated  or  predicted  in  the  simulation  process.  Instead 
of  giving  the  total  current  over  a  certain  area  of  electrolyte 


membrane,  the  average  current  density  is  usually  prescribed 
as  a  basic  given  condition. 

The  activation  polarization  is  a  function  of  the  local 
current  density,  membrane  temperature,  and  gas  pressures, 
which  is  expressed  as 


/(x,y)j?c(x,y) 
io.PoJ^y)  _ 

RTm(x,y )  i(x,y)pll(x,y) 

i  .of  !H2fh2  0>  y ) 


RTm(x,y )  , 

n.Art  =  — — - In 

/act  0.5  F 


(17) 


where  io2  and  /Ho  are  the  constants  of  exchange  current 
densities,  which  are  property  data  of  the  membrane. 

With  the  local  current  density  obtained,  the  local  mass 
fluxes  to  the  fuel  and  air  channels  can  be  correlated  because  the 
electrochemical  reaction  involving  one  mole  of  hydrogen  can 
result  in  proton  flow  (electric  charge)  in  two  moles.  In  the  air 
channel,  oxygen  is  consumed,  and,  as  a  product  of  the  elec¬ 
trochemical  reaction,  water  is  produced.  At  the  same  time,  the 
osmotic  effect  (incorporated  through  the  coefficient  of  ndrag) 
can  cause  a  certain  number  of  water  molecules  to  be  dragged 
from  the  anode  side  of  the  air  channel  through  the  electrolyte 
membrane.  Moreover,  water  might  diffuse  in  between  the 
anode  and  cathode  side  through  the  electrolyte  membrane. 
Considering  these  above  mass  transfer  processes,  the  oxygen 
flux  and  overall  water  flux  on  the  cathode  side  is  formed  as 


mo2  =  Mo2  f 


(18) 
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h2o 


i  i 

^H20  2^  +  47H20^drag  ~ 


nm  _ 

H2°  A7, 
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'm 


(19) 


In  the  fuel  channel,  hydrogen  is  consumed  and  the  balance  of 
water  is  due  to  the  osmotic  effect  and  diffusion  effect: 

<2  =  -MH2  ^  (20) 
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h2o 


IF 

M  „  1  |  r>m  Pity  (+0  ~  +o) 

MH20«drag  p  +  L'fyO  ^  x 
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IT) 


(21) 


where  M  is  the  molecular  weight,  and  2^20  and  are  the 
water  content  of  the  membrane  at  the  interfaces  of  the  anode 
and  cathode  sides.  According  to  Springer  et  al.  [9],  water 
content  can  be  calculated  from  the  activity  of  water  by  the 
correlation: 


J  0.043  +  17.81a  -  39.85a2  +  36.0a3, 
l  14.0+  1.4(a-  1.0), 


where  a  is  given  as 
Xh.ofO,?) 


a  = 


^Sat 

P  h20 


in  which  X  is  the  mole  fraction. 


forO  <  a  <  1 
for  1  <  a  <  3 
(22) 

(23) 
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In  this  model,  we  did  not  consider  the  possibility  of  large 
excesses  of  liquid  water  in  the  flow  channels  or  within  the 
electrodes.  For  the  continuity  equations,  no  distinction  was 
made  between  water  gas  and  liquid  molecules.  It  is  assumed 
that  the  vapor  and  membrane  phases  are  in  phase  equili¬ 
brium.  The  activity  of  water  in  the  vapor  phase  on  anode 
side  is  set  equal  to  the  activity  of  water  in  the  membrane 
[19]. 

The  membrane  properties  of  proton  conductivity,  crm, 
water  diffusivity,  and  osmotic  coefficient,  ndmg,  all  depend 
on  the  water  content  of  the  membrane.  Empirical  expres¬ 
sions  for  the  Nafion  membrane  were  recommended  in 
different  sources  [20];  however,  they  are  mostly  based  on 
the  correlations  by  Springer  et  al.  [9]: 

D'lo  =  Dx  exP  (24 1 6  j  (24) 


where 

'  lO"10, 

D  =  I  10-10(1  +  2(/l  —  2)), 

A  )  10-10(3- 1.67(2 -3)), 
,  1.25  x  IQ-10, 


for  2  <  2 
for  2  <  2  <  3 
for  3  <  2  <  4.5 
for  2  >  4.5 


(25) 


(0.005142 


0.00326)  exp  (^1268 


ndrag  =  0.002922  +  0.052  -  (3.4  x  10“ 19) 


x  102 


(26) 

(27) 


It  may  be  a  disadvantage  that  there  has  been  no  application 
of  these  equations  for  a  temperature  level  lower  than  30  °C. 

The  average  current  density,  initial  gas  species  and  mass 
fractions  are  all  specified  as  known  conditions  of  the  numer¬ 
ical  simulation.  The  fuel  and  air  flow  rates  fed  to  the  fuel  cell 
can  be  determined  through  the  following  two  equations: 


lA  E>Tfuel 
x-'rfuel  m  Aiin 
^in  - 


2  F  pfudXHl 


in  in 


Gair  =  iA  RT?« 


in 


4 F  pair  V^2 
1  in  ^in 


(28) 


(29) 


Table  1 

Fuel  and  air  conditions11 


Mole  fraction 

Air  (1.1  x  105  Pa) 

Oxygen 

0.18 

Nitrogen 

0.81 

Water  vapor 

0.01 

Fuel  (1.1  x  105  Pa) 

Hydrogen 

0.80b 

0.60c 

Water  vapor 

0.20b 

0.40c 

a  Temperature:  85  °C. 
b  Case  1. 
c  Case  2. 


Fig.  3.  Flow  chart  of  the  computation  procedure. 


Fig.  4.  Cell  voltage  vs.  current  density. 
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Potential  Difference  (V) 


X(m) 


0.9069 
0.9031 
0.8993 
0.8955 
0.8917 
0.8879 
0.8841 
0.8803 
0.8764 
0  8726 
0.8688 
0.8650 


Fig.  5.  Local  distribution  of  potential  difference  (E  —  rjact). 


Other  related  parameters  and  property  data  for  the  Nation 
electrolyte  membrane  to  be  used  in  this  work  are  available  in 
the  literature  [12].  These  are: 

•  thickness  of  dry  membrane  (Nafion  112:2mil)  (pm):  51; 

•  density  of  dry  membrane  (g/cm  ):  1.84; 

•  equivalent  weight  of  membrane  (g/mol):  1100; 

•  exchange  current  density  (A/m  ):  /H2  =  100.0  (anode), 
io2  =  1000.0  (cathode); 

•  permeability  of  porous  electrode  (m2):  /?  =  2.3  x  109; 

•  porosity  of  porous  electrode:  s  =  0.4. 

s 
Sr 


2.2.  Numerical  treatment  and  computation  procedures 

A  non-uniform  mesh  arrangement  is  deployed  in  the 
computation  domain,  in  which  the  gas  channels  and  mem¬ 
brane  regions  are  allocated  a  finer  mesh.  The  SIMPLE 
algorithm  [14]  is  adopted  to  link  the  velocity  and  pressure. 
Flow,  temperature  and  mass  fractions  for  each  species  are 
solved  in  each  iteration  step. 

Considering  the  situation  of  one  fuel  cell  unit  in  a  PEM 
fuel  cell  stack,  its  anode  plate  contacts  the  cathode  plate 
of  the  neighboring  unit  while  its  cathode  plate  contacts  the 


Current  Dens  ty  (A/m  ) 


26E-03 
58E  03 
90E-03 
22E-03 
55E  03 
87E-03 
19E-03 
51E-03 
83  E  -  03 
16E-03 
4SE-03 
80E-03 


Xfm) 


Fig.  6.  Local  distribution  of  current  density. 
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anode  plate  of  the  other-side  neighboring  unit.  Suppose  that 
every  fuel  cell  unit  operates  at  the  same  conditions.  Then  it 
is  plausible  to  assume  that  the  outside  surface  temperature 
of  cathode  and  anode  plates  must  comply  to  the  condition  of 
a  consecutive  heat  flux,  that  is: 

q(x,  y )  |z=0  =  q(x,  y) \Z=H  (30) 


The  computational  outflow  boundary  conditions  of  tempera¬ 
ture  and  velocities  for  both  flows  of  fuel  and  air  are  treated 
based  on  the  methods  from  [21,22]. 

The  iterative  computation  follows  the  flow  chart  as  shown 
in  Fig.  3,  and  can  be  described  concisely  as  follows: 

(1)  The  latest  available  mass  fractions  and  pressures  are 
used  to  calculate  the  local  current  densities  and 


Fig.  7.  Velocity  vectors  in  the  air  and  fuel  channels  at  the  half  channel  height:  (a)  fuel;  (b)  air. 
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therefore  the  related  heat  generation,  as  well  as  the 
mass  fluxes  perpendicular  to  the  flow  channels  at  the 
membrane  surfaces.  An  inner  iteration  in  this  step  has 
to  be  conducted  to  guess  and  find  a  proper  cell  voltage, 
so  that  the  integrated  total  current  from  the  local 
current  density  satisfactorily  agrees  with  the  prescribed 
total  current. 

(2)  Equations  of  velocity,  temperature  and  mass  fractions 
are  solved  in  series. 

(3)  The  convergence  of  both  electrical  results  and  the 
results  of  velocity,  temperature  and  mass  fractions  are 


checked.  In  case  no  convergence  is  attained,  additional 
computational  iteration  will  be  necessary,  thereby 
returning  to  the  Step  (1). 

3.  Application  computation 

The  investigated  fuel  cell  unit  has  dimensions  marked  in 
Fig.  1.  Operating  conditions  are  listed  in  Table  1.  With  these 
prescribed  temperatures  and  species’  mass  fractions  in  fuel 
and  air,  the  fuel  cell  electrical  performance  and  local  details 


Mass  Fraction  of  Water 
in  Fuel  Channel 


? 0.0175 


0.015- 


0.0125 - 


0.01- 


0.0075 


I 


0.8551 

0.8364 

0.8178 

0.7992 

0.7805 

0.7619 

0.7432 

0.7246 

0.7059 

0.6873 

0.6686 

0.6500 


0.0075  0.01 

(a) 


0.0125  0.015  0.0175 

X(m) 


Mass  Fraction  ot  Hydrogen 


0.0075 
0.0075 


0.0125 


0.015 


(b) 


0.0175 

X(m) 


0.3171 

0.3010 

0.2849 

0.26S8 

0.2527 

0.2366 

0.2205 

0.2044 

0.1883 

0.1722 

0.1561 

0.1400 


Fig.  8.  Mass  fraction  variation  of  gases  in  the  fuel  channel:  (a)  water;  (b)  hydrogen. 
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of  the  velocity,  mass  fraction,  and  temperature  are  obtained 
and  will  be  presented  in  this  section. 

3.1.  Electrical  performance 

Fig.  4  shows  the  fuel  cell  terminal  voltage  against  current 
density  obtained  through  numerical  calculations  for  series  of 
cell  current  densities.  The  V-I  data  of  a  single  cell  from  the 
work  of  Springer  et  al.  [9]  in  the  similar  membrane  condi¬ 
tions  to  the  present  work  are  plotted  for  reference,  of  which 


the  operation  temperature  was  80  °C  and  Nation  117  mem¬ 
brane  thickness  was  50  pm.  The  cell  voltages  were  found  in 
the  level  quite  approximate.  However,  the  present  fuel  and 
air  channel  structure  is  quite  different  from  that  of  Springer 
et  al.  [9].  This  can  be  the  reason  that  in  the  V-I  curves  there 
are  some  differences. 

Prediction  computation  was  done  for  two  different  cases 
in  different  fuel  quality  or  humidity.  A  significant  effect  of 
humidity  to  the  cell  voltage  can  be  seen  when  current  density 
becomes  larger  as  illustrated  by  the  predicted  curves  of 


Mass  Fraction  of  Oxygen 


^ 0.0175 


0.015- 


0.0125- 


0.01- 


0.0075 
0.0075 


(a) 


0.0125 


0.015 


7) 
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0.0175 
X(m ) 


0.1966 

0.1924 

0.1881 

0.1839 

0.1797 

0.1754 

0.1712 

0.1669 

0.1627 

0.1585 

0.1542 

0.1500 


Mass  Fraction  of  Water  Vapor 
in  Air  Channel 


0.0870 

0.0796 

0.0721 

0.0647 

0.0572 

0.0497 

0.0423 

0.0348 

0.0274 

0.0199 

0.0125 

0.0050 


Fig.  9.  Distribution  of  the  gas  species’  mass  fractions  in  the  air  channel:  (a)  oxygen;  (b)  water;  (c)  nitrogen. 
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Fig.  9.  ( Continued ). 


voltage  versus  current  density.  When  the  mole  fraction  of 
water  vapor  increases  from  20  to  40%,  the  cell  voltage 
increment  can  be  as  large  as  0. 1  V  at  the  current  density  of 
500  mA/cm2. 

The  local  distribution  of  cell  potentials  (E  —  ^act)  at  a 
current  density  of  400  mA/cm  and  fuel  composition  in 
Case  2  is  given  in  Fig.  5.  A  high  potential  can  be  found 
at  the  upstream  region  of  fuel  and  air.  It  then  decreases  along 
the  stream,  and  the  decrement  at  the  upstream  region  is  faster 
than  that  at  downstream  region.  The  corresponding  local 
distribution  of  current  density  is  given  in  Fig.  6.  A  dramatic 
variation  or  decrease  of  the  local  current  density  along  the 
stream  could  be  observed.  At  the  end  of  the  flow  stream  the 
electrochemical  reaction  is  quite  weak. 

3.2.  Velocity  fields 

The  velocity  vectors  in  v-y-plane  at  the  half  height  of  air 
and  fuel  channels  are  shown  in  Fig.  7  for  the  condition  of  a 
current  density  of  400  mA/cm"  and  fuel  composition  of  Case 
2.  In  the  fuel  gas  channel,  velocities  decrease  along  the 
streamwise  direction.  This  reflects  the  depletion  of  mass  in 
the  fuel  channel,  which  includes  the  consumption  of  hydro¬ 
gen  and  loss  of  water  due  to  the  osmotic  effect  that  drags 
water  from  the  anode  to  the  cathode  in  the  electrochemical 
reaction  process.  Although  there  is  water  diffusion  from  the 
cathode  side  to  the  anode  side,  it  is  relatively  weak  and 
cannot  compensate  for  the  water  loss  from  the  osmotic 
effect. 

In  the  air  channel,  consumption  of  oxygen  accompanies 
the  production  of  water  in  the  electrochemical  reaction. 
Contrast  to  the  fuel/anode  channel,  the  air/cathode  channel 
has  a  net  income  of  water  that  is  induced  collectively  from 


the  electrochemical  reaction,  osmotic  drag  effect  and  diffu¬ 
sion  effect.  However,  compared  to  the  overall  mole  flow  rate 
(including  oxygen,  nitrogen  and  water)  in  the  air  channel, 
this  water  is  still  only  a  small  proportion.  This  is  the  reason 
that  we  could  not  see  any  significant  increase  in  the  velocity 
in  air  channel. 

3.3.  Fields  of  mass  fraction 

The  local  distributions  of  the  species’  mass  fractions  in 
the  fuel  channel  for  the  forgoing  discussed  case  are  shown  in 
Fig.  8.  It  is  known  that  both  hydrogen  and  water  vapor  in  the 
fuel  channel  are  consumed  or  transferred  to  the  cathode  side. 
As  a  collective  result  of  the  osmotic  effect  and  mass  diffu¬ 
sion  through  the  membrane,  the  mass  fraction  of  water  vapor 
in  the  fuel  channel  decreases  along  the  streamwise  direction. 
Interestingly,  it  can  also  be  seen  that  the  hydrogen  mass 
fraction  has  a  relative  increment  along  the  stream.  Examin¬ 
ing  the  osmotic  drag  coefficient,  one  can  find  that  its  value 
may  range  from  0.2  to  0.9  depending  on  the  local  water 
content  of  membrane;  therefore,  the  electrochemical  reac¬ 
tion  using  1  mol  of  hydrogen  might  cause  0.4-1. 8  mol  of 
water  to  be  dragged  to  the  cathode  side.  This  is  a  dramatic 
decrease  of  the  mass  fraction  of  water  in  the  fuel  channel, 
and  thus  the  mass  fraction  of  hydrogen  relatively  increases, 
although  the  hydrogen  has  a  net  reduction  due  to  the 
electrochemical  reaction.  From  this  observation,  it  is 
believed  that  water  loss  on  the  fuel  side  is  serious  and  the 
probability  of  water  condensation  in  the  fuel  channel  is  very 
low. 

In  the  air  channel,  oxygen  is  consumed  and  water  is 
generated  due  to  the  electrochemical  reaction  and  osmotic 
effect.  These  variations  can  be  observed  in  Fig.  9,  which 
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displays  the  local  distributions  of  the  mass  fractions  for  each 
of  the  species  along  the  streamwise  direction.  Oxygen  is 
consumed,  and  thus  its  mass  fraction  decreases  along  the 
streamwise  direction.  As  a  product  of  the  electrochemical 
reaction,  water  is  generated  in  the  air  channel.  In  addition, 
the  osmotic  effect  in  the  electrochemical  reaction  drags 
water  from  the  anode  to  cathode,  although  water  inversely 
diffuses  from  the  cathode  to  anode  side.  However,  the  net 
result  is  that  water  increases  in  the  air  channel,  which  can  be 
seen  in  Fig.  9b.  Because  the  nitrogen  and  oxygen  are  in 
rather  high  proportions,  the  partial  pressure  of  the  water 
vapor  is  relatively  low,  and  thus  there  is  no  condensation  of 
water  under  the  investigated  temperature  and  concentration 
conditions. 

Fig.  9c  is  the  mass  fraction  variation  of  nitrogen.  Although 
the  quantity  of  nitrogen  in  air  has  no  variation  in  concentra¬ 
tion,  its  mass  fraction  decreases  because  water  is  produced  in 
air  channel. 

A  careful  examination  of  Figs.  8  and  9  may  find  that  the 
mass  diffusion  in  the  fuel  and  air  channels  exhibits  some 
differences.  In  the  fuel  channel,  the  contour  line  is  rather  flat 
normal  to  flow  direction,  which  indicates  that  the  diffusion 
of  hydrogen  and  water  vapor  is  rather  strong,  and  the  mass 
fraction  is  relatively  uniform  in  the  cross  section  normal  to 
the  flow  stream.  However,  in  the  air  channel,  the  contour 
profile  is  mostly  in  a  parabolic  shape.  This  is  a  reflection  that 
the  diffusion  of  oxygen,  nitrogen  and  water  vapor  is  weak 
compared  to  that  of  the  gases  in  the  fuel  channel,  which  also 
indicates  the  three-dimensional  effect  of  the  transportation 
problem  for  the  cathode  side. 

4.  Concluding  remarks 

A  three-dimensional  numerical  model  simulating  the 
performance  of  a  proton  exchange  membrane  (PEM)  fuel 
cell  with  complicated  gas  channels  is  developed  to  show  the 
strong  coupling  of  the  flow  and  heat/mass  transfer  with  the 
electrochemical  reaction.  The  numerical  model  is  applied  to 
simulate  a  single  PEM  fuel  cell  unit.  The  fields  of  flow,  mass 
fraction,  and  local  electromotive  force  were  obtained  and 
discussed.  The  simulation  results  show  that  the  fuel  gas 
humidity  has  a  significant  effect  on  the  output  voltage  of  the 
fuel  cell.  It  was  also  found  that  due  to  the  electrochemical 
reaction  and  osmotic  effects,  the  air  channel  has  an  increasing 


mass  fraction  of  water,  and  therefore,  the  oxygen  mass 
fraction  or  partial  pressure  is  decreased.  While  high  humid¬ 
ity  of  the  fuel  gas  is  desirable,  excessive  increase  of  water 
in  the  air  channel  could  reduce  the  cell  performance.  An 
optimized  condition  of  humidity  might  exist,  however,  but 
heavy  computations  are  necessary.  It  was  also  found  that 
mass  diffusion  is  strong  in  the  fuel  channel  but  is  rather  weak 
in  the  air  channel.  To  facilitate  a  better  diffusion  of  oxygen 
in  the  air  channel,  modification  of  the  air  channel  and  flow 
arrangement  might  be  ways  to  improve  the  PEM  fuel  cell 
performance. 
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